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Abstract 

By using theoretical analysis and molecular dynamics simulations, we investigate the structure of col- 
loidal crystals formed by nonmagnetic microparticles (or magnetic holes) suspended in ferrofluids (called 
inverse ferrofluids), by taking into account the effect of polydispersity in size of the nonmagnetic micropar- 
ticles. Such polydispersity often exists in real situations. We obtain an analytical expression for the inter- 
action energy of monodisperse, bidisperse, and polydisperse inverse ferrofluids. Body-centered tetragonal 
{bet) lattices are shown to possess the lowest energy when compared with other sorts of lattices, and thus 
serve as the ground state of the systems. Also, the effect of microparticle size distributions (namely, poly- 
dispersity in size) plays an important role in the formation of various kinds of structural configurations. 
Thus, it seems possible to fabricate colloidal crystals by choosing appropriate polydispersity in size. 
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I. INTRODUCTION 



In recent years, inverse ferrofluids with nonmagnetic colloidal microparticles suspended in a 
host ferrofluid (also called magnetic fluid)^"^ have drawn considerable attention for its potential 
application in its industrial applications and potential use in biomedicine.^"^^ The size of the non- 
magnetic microparticles are about 1 ~ 100 fxm, which can be easily made in experiments, such as 
polystyrene microparticles. The inverse ferrofluid system can be modelled in a dipolar interaction 
approximation. Here, the dipolar interaction approximation is actually the first-order approxima- 
tion of multipolar interaction. Because the nonmagnetic microparticles are much larger than the 
ferromagnetic nanoparticles in a host ferrofluid, the host can theoretically be treated as a uniform 
continuum background in which the much larger nonmagnetic microparticles are embedded. If 
an external magnetic field is applied to the inverse ferrofluid, the nonmagnetic microparticles sus- 
pended in the host ferrofluid can be seen to posses an effective magnetic moment but opposite in 
direction to the magnetization of the host ferrofluid. As the external magnetic field increases, the 
nonmagnetic microparticles aggregate and form chains parallel with the applied magnetic field. 
These chains finally aggregate to a column-like structure, completing a phase transition process, 
which is similar to the cases of electrorheological fluids and magnetorheological fluids under exter- 
nal electric or magnetic fields. The columns can behave as different structures like body-centered 
tetragonal (6ct) lattices, face centered cubic (/cc) lattices, hexagonal close packed (hep) lattices, 
and so on. In this work, we assume that the external magnetic field is large enough to form differ- 
ent lattice structures. The actual value of the external magnetic field needed to form such structures 
is related to the volume fraction of the nonmagnetic microparticles and the magnetic properties of 
the host ferrofluid. 

In this work, we shall use the dipole-multipole interaction model^^ to investigate the structure 
of inverse ferrofluids. In ref 12, Zhang and Widom discussed how the geometry of elongated 
microparticles will affect the interaction between two droplets, and introduced higher multipole 
moments' contribution by using a dipole and multipole (dipole-multipole) interaction model to 
give a more exact expression of interaction energy than using the dipole and dipole (dipole-dipole) 
interaction model. The leading dipole-dipole force does not reflect the geometry relation between 
the microparticles nearby, while the dipole-multipole model includes the contributions from the 
size mismatch and is simpler and practical than the multipolar expansion theory ^'^'^'^ in dealing with 
the complex interaction between microparticles for its accuracy. Size distributions can be regarded 
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as a crucial factor which causes depletion forces in colloidal droplets. ^ Even though researchers 
have tried their efforts to fabricate monodisperse systems for obtaining optimal physical or chemi- 
cal properties/^'^^ polydispersity in size of microparticles often exists in real situations/^^^^ since 
the microparticles always possess a Gaussian or log-normal distribution. Here we consider size 
distributions as an extra factor affecting the interaction energy. Polydisperse ferrofluid models are 
usually treated in a global perspective using chemical potential or free energy methods,®'^^'^"^ while 
the current model concerns the local nature in the crystal background. A brief modelling is carried 
out for the size distribution picture in the formation of crystal lattices. The purpose of this paper 
is to use this model to treat the structure formation in monodisperse, bidisperse, and polydisperse 
inverse ferrofluids, thus yielding theoretical predictions for the ground state for the systems with 
or without microparticle size distributions (or polydispersity in size). It is found that when the 
size mismatch is considered between the microparticles, the interaction between them becomes 
complex and sensitive to the different configurations. This method can also be extended to other 
ordered configurations in polydisperse crystal systems. 

This paper is organized as follows. In Section UIl based on the dipole-multipole interaction 
model, we present the basic two-microparticle interaction model to derive interaction potentials. 
In Sections [Inland |IVl we apply the model to three typical structures of colloid crystals formed in 
inverse ferrofluids, and then investigate the ground state in two different configurations by taking 
into account the effect of size distributions. As an illustration, in Section |V] we perform molecular 
dynamics simulations to give a picture of the microparticle size distribution in the formation of a 
bet lattice in bidisperse inverse ferrofluids. The paper ends with a discussion and conclusion in 
Section IVIl 

II. INTERACTION MODEL FOR TWO NONMAGNETIC MICROPARTICLES 

We start by considering a simple situation in which two nonmagnetic spherical microparticles 
(or called magnetic holes) are put nearby inside a ferrofluid which is homogeneous at the scale of 
a sphere in an applied uniform magnetic field H, see Figure 1 . The nonmagnetic microparticles 
create holes in the ferrofluid, and corresponding to the amount and susceptibility of the ferrofluid, 
they possess the effective magnetic moment, which can be described by^^ m = —j^r^VM = 
—xVH, where Xf (or X) means the magnetic susceptibility of the host (or inverse) ferrofluid. 
When the two nonmagnetic microparticles placed together with distance away, we can view 



3 



the magnetization in one sphere (labelled as A) is induced by the second (B). The central point 
of dipole-multipole technique is to treat B as the dipole moment m at the first place, and then 
examine the surface charge density S induced on the sphere A. From S we can use the multipole 
expansion (detailed discussion can be found in ref 26) to obtain the multipole moment. When 
exchanging the status of A and B, treating A as the dipole moment, the averaged force between 
the two microparticles is thus obtained. For the perturbation of the magnetic field due to the 
two microparticles, magnetization M in the microparticles become nonuniform, and they will 
obtain multipole moments from mutual induction. However, the bulk magnetic charge density still 
satisfies p = V ■ M = 0. So we need only study the surface charge S^i, 

S„ = n ■ M, (1) 

where fi is the unit normal vector pointing outwards. The magnetic multipole moments by surface 
charge density in spherical coordinates can be written as, 

qin = j l^:(0,(/p)rlS„d5, (2) 

where ri denotes the radius of the microparticle. All n 7^ moments vanish due to rotational 
symmetry about the direction of magnetization, so eq 2 can be rewritten as 

Ai = qio = V(2; + l)7r j P;(cos(/))rl+2S„(cos0)dcos0. (3) 

We can expand the surface charge density in Legendre polynomials in the spherical coordinates 

(r,0,v9),^2 

Mcos<P). (4) 



00 X \ I— I 

S„(cos0) = -J^(_l)m;(^ + 1) 11 p^, 



The / > 2 parts in eq |4] correspond to the effects of multipole (that are beyond the dipole). Here 
we set spherical harmonics dy? ^ sin 0^*^/ (0, V^)d0 = 5^' „. The force between 

the dipole moment m and induced multipole moment Ai can be derived as 

Fd-m = + + 2)(^)^^ COS0. (5) 

In view of the orthogonal relation Jj^^ PpPpdx = 2!^! (!-")! "^^'^ obtain the interaction energy 
for the dipole-multipole moment 

7T TT /(/ + 1)2 2/+1 1 ~ 3 cos^ 6' l-3cos^6 
Ud-m = ^/i/mim2 }^ -X \^ [ {rf+^ + rf+^) 2lT^ = /^Z^ /(O 2lT^' (6) 

I ^ij I ^ij 



where // = lnfmim2, f{l) = f X^§q^ + r^'+^), the suffix D (M) of force F or energy U 
stands for the dipole moment (multipole moment), and magnetic permeabiHty /i/ = //o(l + X/) 
with /io = Att X 10~^ H-m^^ Here mi and m2 denote the effective magnetic moments of the two 
nonmagnetic microparticles, which is induced by external field as dipolar perturbation. < 6' < | 
is the angle between the their joint line with the direction of external field, and and (p are both 
the spherical coordinates (r,(f),(p) for one single nonmagnetic microparticle. For typical ferrofluids, 
there are magnetic susceptibilities, Xf — 1-9 X — 0.836.^^ Because we consider the bet, fee 
and hep lattices, the crystal rotational symmetry in the xy plane is fourfold, and the value of n 
can only be 0, ±4, .... In the general case, when the polarizabilities between the microparticles 
and ambient fluid is low, the higher magnetic moments can be neglected since they contribute less 
than 5 percent of the total energy. In this picture, the nonmagnetic microparticle pair reflects 
the dipole-dipole and dipole-multipole nature of the interaction, and can be used to predict the 
behavior of microparticle chains in simple crystals. 

III. POSSIBLE GROUND STATE FOR UNIFORM ORDERED CONFIGURATIONS 

Let us first consider a bidisperse model which has been widely used in the study of magnetorhe- 
ological fluids and ferrofluids. The model has large amount of spherical nonmagnetic micropar- 
ticles with two different sizes suspended in a ferrofluid which is confined between two infinite 
parallel nonmagnetic plates with positions at z — and z — L, respectively. When a magnetic 
field is applied, dipole and multipole moments will be induced to appear in the spheres. The in- 
verse ferrofluid systems consist of spherical nonmagnetic microparticles in a carrier ferrofluid, and 
the viscosity of the whole inverse ferrofluid increases dramatically in the presence of an applied 
magnetic field. If the magnetic field exceeds a critical value, the system turns into a solid whose 
yield stress increases as the exerting field is further strengthened. The induced solid structure is 
supposed to be the configuration minimizing the interaction energy, and here we assume first that 
the microparticles with two different size have a fixed distribution as discussed below. 

Using the cylindrical coordinates, the interaction energy between two microparticles labelled 
as i and j considering both the dipole-dipole and dipole-multipole effects can be written as 

rr / N /-.V- /(O N ,1 -3C0S2^, 

Uijip, Z) = //(I + J] 4^) • ( -3 ), (7) 

I "^ij »J 

where the center-to-center separation rij — | — r ^ | = [p^-\-{zi — ZjY]^, and ^ is the angle between 
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the field and separation vector rij (see Figure 1). Here p — [{xi — XjY + (t/j — yjf] ^ stands for the 
distance between chain A and chain B (Figure 2), and Zi denotes the vertical shift of the position 
of microparticles. Since the inverse ferrofluid is confined between two plates, the microparticle 
dipole at [x, y, z) and its images at (.x, y, 2Lj±z) for j = ±1, ±2, ... constitute an infinite chain. In 
this work, we would discuss the physical infinite chains. After applying a strong magnetic field, the 
mismatch between the spheres and the host ferrofluid, as well as the different sizes of the two sorts 
of spheres will make the spheres aggregate into lattices like a bet (body-centered tetragonal) lattice. 
In fact, the bet lattice can be regarded as a compound of chains of A and B, where chains B are 
obtained from chains A by shifting a distance ri (microparticle radius) in the field direction. Thus, 
we shall study the case in which the identical nonmagnetic microparticles gather together to form 
a uniform chain, when phase separation or transition happens. For long range interactions, the 
individual colloidal microparticles can be made nontouching when they are charged and stabilized 
by electric or magnetic static forces, with a low volume fraction of nonmagnetic microparticles. 
The interaction energy between the nonmagnetic microparticles can be divided into two parts: 
one is from the self energy of one cham(Us), the other is from the interaction between different 
chains(f/jj(p, z)). Consider the nonmagnetic microparticles along one chain at rj = 2ajz {j = 
0, ±1, ±2, ...) (namely, chain A), and the other chain at rj = {2j + l)az (chain B), the average 
self energy per microparticle in an infinite chain is Ug — —n S^i[ (2as)3 + 2 (2^p+4]- 

If we notice that for an infinite chain all even mulitipole contributions vanish due to spatial 
magnetic antisymmetry around the spheres, the sum starts at Z = 3. Because the radius of the 

^2i+l , 2i+l 

sphere is smaller than the lattice parameter o, for large multipole moment, ^ ^^^-^i — ^ 1, we 

need only consider the first two moment contributions for simplicity. Thus the average self energy 

[/. can be calculated as = -2^(0 + «P + Z^) = -^(0« + omn^jom + 

0.00195507/(5) -J ^ whcrc ({n) — ^ is Riemann ( function. The interaction energy between 

two parallel infinite chains can be given by \Uij{p, z), in which the microparticles along one chain 

locate at Vj — 2ajz {j — 0, ±1, ±2, ...) and one microparticle locates at rj — p + zz, 

9 °° 1 
U,j{p,z) = -^[(2 + p— ) J2 



9p jt'o.ip' + i^-'^j^m 

3 d °° 1 

; j=—oo '■ ^ ^ / J 

U, + U2. (8) 



Following the Fourier expanding technique which is proposed by Tao et al.,^'^ we derive U2 




+ 2-2-^{^y+ln'+'cos{ — )-S] 
a a 



(9) 



with 



s=l 
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Here Ki{x) represents the ith order modified Bessel function, r(x) the T function, and s denotes 
the index in Fourier transformation.^^ And the dipole-dipole energy Ui is written as 



We obtain the expression for Uij{p, z), and the interaction energy per nonmagnetic microparti- 
cle U{p, z) is [/s + I Yl,k ^ijiP-' where Yl,k denotes the summation over all chains except the 
considered microparticle. For the same reason of approximation discussed above, we need only 
choose the first two terms (/ = 3 and / = 5) in the calculation. 

The interaction between chain A and chain B depends on the shift z, the lattice structure 
and the nonmagnetic microparticle size. An estimation of the interaction energy per nonmag- 
netic microparticle includes the nearest and next-nearest neighboring chains, here we could 
discuss three most common lattice structures: hct, fee, and hep lattices. For the above lat- 
tices, their corresponding energy of Uij{p, z) can be respectively approximated as Uij^bctip, z) = 
Wij{V^a, ^ = 0) - 4[/,j(v^a, z = 0), U.jjcdp, z) = 4Uij{V3a, z = 0) - 2Uij{2a, z = 0), and 
Uij^hcpip, z) = 3Uij{^/3a, z = 0) - AUij{2a, ^ = 0). 

Figure |3] shows, for different lattices, the dependence of Uij{p, z) on the vertical position shift 
z, which determines whether the interaction is attractive or repulsive. Uij{p, z) reflects the energy 
difference between chain A and chain B for (a) bet, (b) fee, and (c) hep lattices. It is evident 
that, for the same lattice structure, t/jj(p, z) is minimized when the size difference between chain 
A and chain B is the smallest. For the sake of comparison, we also plot the results obtained by 
considering the dipole-dipole interaction only. Comparing the different lattices, we find that the 
hct lattice possesses the smallest energy at the equilibrium point, thus being the most stable. 

Figure m displays the interaction energy U{p, z) as a function of the lattice constant a for the 
hct lattice. It is shown that as the lattice constant increases, the dipole-multipole effect becomes 
weaker and weaker, and eventually it reduces to the dipole-dipole effect. In other words, as the 
lattice constant is smaller, one should take into account the dipole-multipole effect. In this case, 
the effect of polydispersity in size can also play an important role. 



oo 
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Figure \5\ displays the interaction energy per nonmagnetic microparticle U {p, z) vs the lattice 
parameter a for different lattice structures. The bet structure also proves to be the most stable state 
while the hep lattice has the highest energy. It also shows that the energy gap between hot lattice 
and fee lattice exists but is small. Figure|6ta) shows that the energy gap At/ = U^ct — Ujcc is about 
0.5 percent of the interaction energy value. In this aspect, the bet lattice proves always to be a more 
stable structure comparing with fee. As the radius of microparticles increases, the energy gap 
between bet and fee lattice enlarges accordingly. That is, the bet lattice becomes much more stable. 
Figure |6tb) shows the bet lattice energy U {p, z) in respect of different sizes of microparticles 
for chain A and chain B. It can be seen that the close touching packing (ri = r2 = a) has 
the lowest energy state. However, also from the graph, the crystal with the same microparticle 
size (monodisperse system) may not be the lowest energy state, which gives a possible way of 
fabricating different crystals by tuning the distribution of microparticle size. 

IV. POLYDISPERSE SYSTEM WITH RANDOM DISTRIBUTIONS 

In Section Unl we have discussed the structure and interaction in a bidisperse inverse ferrofluid 
(namely, containing microparticles with two different sizes). But the interaction form in poly- 
disperse crystal system is complex and sensitive to the micro structure in the process of crystal 
formation. Now we investigate the structure of polydisperse inverse ferrofluids with microparti- 
cles of different sizes in a random configuration. To proceed, we assume that the average radius r 
satisfies the Gaussian distribution 



where o denotes the standard deviation of the distribution of microparticle radius, which describes 
the degree of polydispersity. Integrating eq|9]by ri and r2, we could get the average dipole- 
multipole energy U2- Doing the same calculation to self energy Ug, we can get the average inter- 
action energy U (p, z) = U S + U1 + U2, where the microparticle size ri and r2 are replaced by the 
mean radius tq. The microparticle sizes will be distributed in a wider range as long as a larger a 
is chosen. 

Figure IT^a) shows the ground state interaction energy of bet lattice for different polydisperse 
distributions. As the degree of polydispersity a increases, the energy U (p, z) drops fast, especially 
when the distribution of microparticle size gathers around r = a. It shows that the inverse fer- 
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rofluid crystal in the formation of ground state tends to include microparticles possessing more 
different sizes. The crystal configuration energy of two uniform chains with identical micropar- 
ticle aggregation is also plotted in the graph. Here we consider two cases. First, we assume the 
microparticles in chain A and chain B are identical, ri = r2 = tq. As tq increases, the behavior of 
energy decreasing is discovered to be similar with the random configuration with a = 0.2ro. Sec- 
ond, we set for one chain, such as chain B, the microparticle size r2 = a to be unchanged, while 
the microparticle size ri of chain A increases. It shows that the bet lattice energy for the second 
case is lower than the first case and two other random configurations. And it also shows that the 
random configuration is not always the state with the lowest energy. It proves that polydisperse 
systems are sensitive to many factors which can determine the microstructure. Figure IT^b) shows 
the energy gap between bet and fee lattices for different distribution deviation a. It is evident that 
higher a leads to larger energy difference between bet and fee, especially at larger tq. In other 
words, at larger tq and/or a, bet lattices are much more stable than fee. 



V. MOLECULAR-DYNAMIC SIMULATIONS 

Here we use a molecular-dynamic simulation, which was proposed by Tao et al.,^^ in order 
to briefly discuss the structure formation of bidisperse inverse ferrofluids. The simulation herein 
involves dipolar forces, multipole forces, viscous drag forces and the Brownian force. The mi- 
croparticles are confined in a cell between two parallel magnetic pole plates, and they are ran- 
domly distributed initially, as shown in Figure 8(a). The motion of a microparticle i is described 
by a Langevin equation, 

rria-j-^ = Fi - STvar]—^ + Ri{t), (13) 
at^ at 

where the second term in the right-hand side is the Stokes's drag force, Ri is the Brownian force, 
and 

F^ = Y.(f^j + fin + fr''- (14) 

Here fij = —VU {p, z), while fl^^, f^""'''^ and Ri{t), have the similar expressions as those in ref 27 
and the references therein. In eq 13, ma and a are respectively the average mass and diameter of 
microparticles. Figure [8] shows the inverse ferrofluid structure with the parameters, magnetic field 
H = 14 Oe, temperature T = 300 K, A' = ^^f^m^ = and B' = V^^f'BTa^^/r ^ ^q„4_ 

Here X denotes the ratio of a dipolar force to a viscous force, B' the ratio of Brownian force to a 
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dipolar force, ks the Boltzmann constant, and r the subinterval time step. 

We take into account a bidisperses system that contains two kinds of microparticles with dif- 
ferent sizes, as shown in Figure 8. In details, this figure displays the configuration of microparticle 
distribution in a bidisperse system at (a) the initial state, (b) the state after 15 000 time steps, and 
(c) the state after 80 000 time steps. The order of Figure 8(c) is better than Figure 8(b). Here we 
should remark that 80 000 steps give the sufficient long time steps to reach the equilibrium state 
for the case of our interest. The structure for the bi-disperse system in Figure 8(b) and (c) has 
the following features: (i) In the field direction, the large spheres form the main chains from one 
plate to the other, where the large spheres touch each other, (ii) The large spheres also form many 
small bet lattice grains. However, they do not form a large bet lattice, (iii) The small spheres fill 
the gaps between these bet lattice grains. From Figure 8, it is observed that, for the parameters 
currently used, the order of a bidisperse system (which is a bct-]ike structure) is not as good as that 
of monodisperse system (no configurations shown herein). We should remark that the long-range 
interaction can yield the above-mentioned bet lattice structure, but some perturbations caused by 
the Brownian movement existing in the system can change it to another lattice structure which 
has similar free energy. Therefore, for the bi-disperse system of our interest, the large spheres 
form the main chains from one plate to the other in the field direction, thus forming many small 
bct-like lattice structures. While small spheres fill the gaps between these bct-like lattices, they 
themselves do not form a bct-like lattice due to such perturbations. Here we should also mention 
that the degree of order of a specific system depends on the choice of various physical parameters, 
for example, the size of microparticles and so forth. 



VI. SUMMARY 



In summary, by using theoretical analysis and molecular dynamics simulations, we investi- 
gate the structure of colloidal crystals formed by nonmagnetic microparticles (or magnetic holes) 
suspended in a host ferrofluid, by taking into account the effect of polydispersity in size of the non- 
magnetic microparticles. We obtain an analytical expression for the interaction energy of monodis- 
perse, bidisperse, and polydisperse inverse ferrofluids. The bet lattices are shown to possess the 
lowest energy when compared with other sorts of lattices, and thus serve as the ground state of the 
systems. Also, the effect of microparticle size distributions (namely, polydispersity in size) plays 
an important role in the formation of various kinds of structural configurations. Thus, it seems 
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possible to fabricate colloidal crystals by choosing appropriate polydispersity in size. As a matter 
of fact, it is straightforward to extend the present model to more ordered periodic systems,^^ in 
which the commensurate spacings can be chosen as equal or different. 
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Figure captions 

Figure 1. Schematic graph showing two nonmagnetic microparticles (magnetic hole) of radius 
ri and r2, suspended in a ferrofluid under an applied magnetic field H. 

Figure 2. Three different lattices, hot, fee, and hep, which are composed of non-touching 
microparticles with different size distribution. 

Figure 3. The dependence of interaction energy z) (in units of /xq) versus vertical shift z 
for different lattices: (a) het, (b) fee, and (c) hep. In the legend, "dipole-dipole" denotes the case 
that the dipole-dipole interaction is only considered for calculating the interaction energy. 

Figure 4. The interaction energy U (p, z) versus lattice constant a. The solid line stands for the 
case in which the dipole-dipole interaction is only considered. 

Figure 5. The interaction energy per microparticle U{p,z) versus lattice constant a for different 
lattices, bet, fee, and hep. 

Figure 6. (a) The energy gap AU — Ubct — Ufcc (in units of /io) versus different size of 
nonmagnetic microparticles in chain A and chain B; (b) The het lattice energy ?7(p, z) (in units of 
Po) versus different sizes of nonmagnetic microparticles in chain A and chain B. 

Figure 7. (a) The ground state interaction energy of a het lattice versus microparticle size Tq 
for random polydispersity configuration (solid, dashed, and dotted lines) and a configuration com- 
posed of two different uniform chains (dash-dotted and short-dash-dotted lines), (b) The energy 
gap At/ between het and fee lattices for random polydispersity versus Tq. 

Figure 8. The configuration of nonmagnetic microparticle distribution at (a) the initial state, 
(b) the state after 15000 time steps, and (c) the state after 80000 time steps. 
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FIG. 1: Jian, Gao, Huang, and Tao 
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FIG. 2: Jian, Gao, Huang, and Tao 
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FIG. 4: Jian, Gao, Huang, and Tao 
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FIG. 6: Jian, Gao, Huang, and Tao 
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FIG. 7: Jian, Gao, Huang, and Tao 



21 




FIG. 8: Jian, Gao, Huang, and Tao 
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